#coding=utf-8
import numpy as np
import matplotlib.pyplot as plt

def MCIntegrate( f, a, b, n ):
    x = a + (b-a) * np.random.random(n)
    samples = f(x)
    fave = np.average(samples)
    f2ave = np.average( samples **2 )
    sigmaprime = np.sqrt( (f2ave - fave**2) / n )
    print( "99.7% possibility: Integral in [",(b-a)*(fave - 3*sigmaprime), ",", (b-a)*(fave + 3*sigmaprime),"]" )
    #print("fave = ", fave, "f2ave = ", f2ave)
    #print("samples = ", samples)

def f(x):
    return np.sin(x)

MCIntegrate( f, 0, np.pi, 100000000)